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INTRODU  CT  ION 


This  is  the  first  of  two  tee* hwc&I  reports  derived  from  the  disserta- 
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tion  by/Yao,  Neng-Chun  titled  ^"Bispectral  and  Cross -bis pectral  Analysis 
of  Wmds  and  Currents  off  the  Oregon  Coast.  " In  this  first  technicaP  PJ  L 

"i  ._v  ^ J5  f O v i (Lf  w— - 

report  we  give^the  theoretical  development  of  the  rotary  component 


approach  to  the  bispectrum  and  cross -bispectrum  among  vector  random 
processes,  as  well  as  the  development  of  a set  of  energy  transfer  func- 
tions by  which  the  nonlinear  interaction  among  triplet  frequency  com- 
ponents of  the  processes  can  be  quantized.  Statistical  tests  are  dis- 


cussed; an  appendix  contains  computational  information. 

voi  j ft  prc  M 'c iff 


In  the  subsequent  technicapreport  we  give  the  results  of  application 
of  the  rotary  component  method  of  cross -bispectrum  analyses  to  a set 
of  winds  and  current  data  acquired  from  the  research  buoy  TOTEM. 
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DEVELOPMENT 


Most  physical  processes  in  the  atmosphere  and  in  the  ocean  are 
not  strictly  Gaussian  in  nature  (Tukey,  1961)  and  the  classical  spectrum 
analysis  does  not  reveal  any  correlation  and  any  energy  transform  due 
to  nonlinear  interactions.  Higher  order  spectrum  analysis  techniques 
are  needed  to  address  such  a problem.  Existing  techniques  of 
bispectrum  and  cross-bispectrum  analysis  are  particularly  useful  for 
studying  the  quadratic  interactions  (lowest  order  of  non-linear  inter- 
actions). However,  they  have  been  developed  along  the  assumption 
that  all  the  processes  are  real.  For  studying  vector  random  processes, 
such  as  wind  and  ocean  current  velocity,  one  is  forced  to  resolve  them 
into  scalar  quantities  along  certain  rectangular  coordinates.  The  results 
of  the  spectrum  analysis  from  the  resolved  data  not  only  depend  highly 
on  the  choice  of  the  coordinates  but  also  are  very  difficult  to  interpret 
(Fofonoff,  1969;  Mooers,  1970). 

A new  version  of  bispectrum  and  cross-bispectrum  analysis, 
based  on  the  rotary  component  concept  given  by  Mooers  (1970),  is 
developed  here  to  deal  specifically  with  this  particular  problem. 


2 







An  energy  transfer  equation  in  terms  of  rotary  spectral  parameters 
is  given.  It  will  show  the  certain  amount  of  energy  transformed 
due  to  particular  quadratic  interactions  involved.  A brief  review  of 
bispectrum  and  cross -bispectrum  analysis  as  well  as  the  rotary 
component  spectrum  analysis  will  be  given  before  the  detail  of  the 
development  of  the  new  technique  is  presented.  The  newly -developed 
technique  will  be  used  to  study  the  quadratic  interactions  between  the 
atmospheric  wind  and  the  near -surface  ocean  current. 


A.  Bispectrum  and  Cross -bispectrum  Analysis 

Since  the  suggestion  made  by  Tukey  (1961),  several  studies 
of  the  estimation  of  higher  order  spectra  have  been  given  (Rosenblatt 
and  Van  Ness,  1965;  Van  Ness,  1966;  Brillinger  and  Rosenblatt, 

1966b;  Hinich  and  Clay,  1968;  and  Nagata,  1970).  Bispectrum 
analysis  has  been  applied  to  the  study  of  non-Gaussian  processes 
in  geophysical  problems  (Hasselmann,  Munk,  McDonald,  1963; 
McDonald,  1963;  Haubrich,  1965;  Murty  and  Henry,  1972)  as  well  as 
in  the  other  fields  (Godfrey,  1965;  Hasselmann,  1966;  Barnett, 
Johnson,  Naitoh,  Hicks  and  Nute,  1971;  Huber  , Kleiner,  Gasser  and 
Dumermuth,  1971).  Application  of  cross -bispect  rum  analysis  has 
been  given  by  Hasselmann's  (1966)  study  of  nonlinear  ship  motion, 
and  by  Roden  and  Bendiner's  (1973)  study  of  oceanic  parameters. 
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Higher  order  spectrum,  such  as  the  trispectrum,  have  not  yet  been 
used  in  practical  problems  because  they  require  large  amounts  of 
computing  time  (Nagata,  1970).  Even  bispectrum  and  cross -bispectrum 
analysis  are  relatively  new.  However,  they  are  useful  in  studying  the 
low-order  nonlinear  interacting  phenomenon  associated  with  non- 
Gaussian  processes. 

A brief  review  of  the  definitions  concerning  their  various  functions 
and  their  physical  meaning  will  be  given  in  the  following  section. 
Comparisons  with  ordinary  spectrum  functions  are  made  for  better 
understanding.  An  effort  is  made  to  emphasize  the  difference 
between  the  spectral  density  function  of  a frequency  and  the  part  of 
spectrum  contained  in  a resolved  frequency  band  around  that  particular 
frequency. * 

1.  Definitions  of  Bispectrum  Density,  Biphase  and  Bicoherence 


The  spectral  and  bispectral  density  of  a random  process  are 

' 


The  physical  meaning  of  spectrum  density  and  the  part  of  the  spectrum 
contained  in  a resolved  frequency  band  centered  at  the  particular  fre- 
quency are  obviously  different.  However,  in  most  literature,  a 
distinction  between  these  two  quantities  has  not  been  established,  and 
the  definition  of  the  spectrum  density  by  Jenkins  and  Watts  (1968)  adds 
further  confusion.  They  define  the  spectrum  density  function  as  the 
normalized  power  spectrum  (i.  e.  , the  part  of  the  spectrum  in  a 
resolved  frequency  band  divided  by  the  total  variance  of  the  random 
process).  For  the  sake  of  consistency,  their  definition  will  not  be 
used  in  this  study. 
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defined  as  the  frequency  decomposition  of  the  random  processes 
second  and  third  moment  respectively  (Hasselmann,  et  al.  , 1963). 

For  example,  if  £(t)  is  a stationary  random  function  of  time  (or  dis- 
tance) with  zero  mean,  the  spectrum  density  F (X)  and  bispectrum  den- 
sity B0y  \2)  of  C (t)  are  obtained  from  the  Fourier  transforms  of 
£(t)'s  mean  second-  and  third-order  products: 


F(M  = ^-  J°°  R(T)  e-i2^Tdr 

- OC 

where  R(T)  = < C(t)C(t  + T)>  , and 


a) 


B'VV=-77  if  S(T1,T2)e-i2"l\Tl+,l2Ti)dT  dT 
(^tt)  -oo 


(2) 


where  S(T,T  ) = < c(t)C(t  +r  K(t  +T  ) >. 

i £ 12 

R(T)  and  S(Tj,  T^)  are  called  the  covariance  and  bicovariance  of  £(t) 
respectively,  j is  the  time  lag  (or  distance  lag),  and  \ is  the  fre- 
quency (wave  number  in  case  that  Tis  the  distance  lag).  The  brackets 
<>  denote  means.  The  inverse  relations  of  (1),  (2)  are 


„ i2-n-X-T,. 

F (\)e  d\ 


R<T»=y 

- 0C 

S(T1-  rz)  J B(\1,\z)ei2lT(KlTl+X2T2,dXid\2  . 


(3) 

(4) 
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Here  * denotes  the  conjugate  of  the  quantity.  As 


(t)  is  assumed  to 


be  stationary  with  real-valued  components,  the  known  symmetry 


relations 


follow  (Has selmann  et  al.  , 1963;  Brillinger  and  Rosenblatt,  1966a,  b) 


In  terms  of  the  spectrum  and  bispectrum  densities,  (7)  and  (8)  become 


From  (5).  (6),  (9  ) and  ( 1 0 ) » it  follows  that  the  spectrum  den 


bispectrum  density  is  complex  (see  below)  and  is  determined  bv 


its  values  in  an  octant;  for  example 


0 < X 2 < oc  . 0 < X.  j < X.  ^ • 


The  spectrum  density  and  bispectrum  density  of  a random  process 
can  also  be  expressed  in  terms  of  the  components  dZ(X)  of  the  Fourier- 
Stieltjes  representation  of  the  random  process  $(t)  as  follows  (Has selm; 
et__al.  , 1 963;  MacDonald,  1963;  Huber _et  al.  , 1971)  ; 

C (t)  = f dZ(\)el2lTKt 


i2irXt 


<dZ(\  )dZ(X  ) > = F(X  )dX,  if  V +X  = 0 
a 1112 


if  X +X2  t 0 , 


i.e.,  < dZ(Xj  )dZ  "'(Xi  ) > = F(\  )dX  , ar 


< dZ(X1)dZ(X2)dZ(X3)>=B(XI,X2)dX1dX2,  if  XJ+  X^  + X = 


if  Xf  + X^  + X3  j-  0 > 


i.e.,  <dZ(X  )dZ(X  )dZ  (X  +X  ) > = B(X  , X )dX  dX 

^ id  12  12 

These  expressions  not  only  provide  an  alternative  computational  pro- 
cedure for  estimating  the  spectrum  and  bispectrum  densities  but  also 
suggest  some  physical  interpretations  of  them.  From  equation  (11), 
one  can  see  that  F (AjJdAj,  the  part  of  spectrum  in  the  resolved 


frequency  band  centered  at  X^,  represents  the  contribution  to  the 

~2 

mean  square  of  the  random  process,  £ (t),  from  the  products  of  its 
Fourier  components  of  frequencies  X^  and  -X^ . From  equation  (12) 
one  finds  that  B(X  , X )dX  dX  the  part  of  bispectrum  in  the  two- 

L 1 ' 

dimensional  frequency  domain  bounded  by  the  resolved  frequency 

bands  centered  at  X^  and  X^  respectively,  represents  the  contribution 

3 

to  the  mean  cube  of  the  random  process,  £ (t),  from  the  product  of 

its  three  Fourier  components  of  frequency  X^ , X^  and  - (X^  +X^)  resPec~ 

tively.  It  is  interesting  to  note  that  bispectrum  density  is  a function 

of  two  independent  frequencies  only.  Since  the  mean  cube,  or  third 

moment,  of  a random  process  generally  results  from  a second  order 

interaction  (Nagata,  1970),  B(X  ,X  ) is  a measure  of  the  amount  of 

1 £ 

such  interaction  involving  the  random  processes  components  whose 
frequencies  are  X^,  X^  and  - (X^  +X^)  respectively.  Since  B(X^ , X^)  is  a 
complex  function,  (as  shown  by  equation  (12)  , the  interpretation  of 
the  interactions  among  the  components  would  be  clearer  if  B(X^.X^) 
is  expressed  in  polar  form. 


B<h>2> 


B<W 


1 lrn[  B(Xj , X 2 ) ] 

+ (W  = tan'  Re[B(K;,r2';) 


(1  3) 


-IT  < <j>(X  , X ) < IT 


_ _ 


| B (X  j , X.  ^ ) I *s  ^le  amplitude  of  the  bispectrum  density.  If  it  is  not 
significantly  larger  than  zero,  there  are  no  interaction?  among  the 
components  of  the  three  frequencies.  4>  (\  x , \ 2)  is  the  phase  of  the 
bispectrum  density.  It  is  determined  by  the  ratio  of  the  imaginary 
part  of  the  bispectrum  density  to  the  real  part  of  it,  and  was  formally 
called  biphase  by  Barnett  et  al.  (1971). 

A third  function  which  is  derived  from  the  normalization  of  the 
bispectrum  has  been  given  various  names  in  literature.  Hasselmann 
et  al.  (1963)  relates  the  normalized  bispectrum  to  the  skewness  with 
the  following  dimensionless  equation 

C2V  B(\  j , \ ? )d\  ^ 

~ 2 2 2 4 = 3 I (14 

[<^  ><?2  ><?3  >]2  [ F(XJ  )F(\2)F(\3)d\'5]2 

where  Cj  • are  the  frequency  components  of  £(t)  obtained  by 

passing  it  through  three  filters  centered  on  V , \ , \ _ and 

of  bandwidth  dk  . Godfrey  (1965)  calls  the  result  of  the  same  equation 
the  bispectral  coefficient.  Since  it  is  a measure  of  skewness  of  the 
random  process,  it  is  also  called  the  bispectral  skewness  (Hinich  and 
Clay,  1968).  The  name  bicoherence  also  appears  in  the  published 
papers  (Haubrich,  1 965;  Bar  nett  et  al. , 1 971 ; Huber  et  al.  , 1971), 
but  it  has  been  defined  rather  inconsistently  by  various  authors. 
Haubrich  (1965)  and  Barnett  et  al.  (1971  ) show 


2 B(X1,X2)  2 

bic  (X  , X ) = r — 

< X(X  X(X  2 ) X(Xj  + X2r  > 


where  bic(X  ^ , X^)  is  the  bi  coherence.  Here  (B(X  ^ , X ^)  is  defined  as 
bispectrum  (rather  than  bispectrum  density)and  B(X  , X = < X(X  ^ ) • 

X (X  „ )X*(X  , + X , ) > and  X(X  . ) is  the  F ourier  coefficient  of  the  random 

2 12  j 

process  obtained  from  discrete  Fourier  transform. 

Equation  (15)  is  essentially  the  same  one  with  which  Hinich  and 
Clay  (1968)  used  to  define  the  bispectral  skewness,  except  they  explicitly 
express  it  as  the  spectrum  over  the  resolved  band  centered  at 
frequency  X.  Equation  (15)  appears  differently  from  equation  (14) 
given  by  Hasselmann  et  al.  (1963)  because  different  scales  of  Fourier 
transform  parameters  are  used.  (See  Hinich  and  Clay,  1968,  for 
details.  ) Physically,  however,  they  are  the  same.  Thus  by  definition, 
the  bicoherence  or  bispectral  skewness  is  a dimensionless  real  value 
with  a range  from  zero  to  one. 

A definition  given  by  Huber  _et  al.  (19  71)  seems  different.  He 
also  calls  the  bicoherence  the  normalized  bispectrum,  but  he  gives  a 
mathematical  relationship  as 


B(X  , X ) 

bic(X  , X ) = r 

(F(X  )F(X  )F(X  ))2 


whe  e B(\  ^ , >>  is  the  bispectrum  density  and  F(  j)  is  the  spectrum 
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density  as  defined  by  equations  (10)  and  (11).  Bicoherence  so 


defined  has  a dimension  of  (dX  ) 2 and  is  not  a real  value.  The 


amplitude  varies  from  zero  to  (dX  ) 2.  For  the  sake  of  consistency. 


this  definition  of  bicoherence  will  not  be  used.  In  this  study,  bico- 


herence will  be  considered  as  the  dimensionless  absolute  value  of  the 


normalized  bispectrum  or  bispectral  skewness  defined  in  equation 


Hinich  and  Clay  (1968)  show  that  if  the  component  at  \ ^ + \ ^ is 


basically  due  to  a multiplicative  interaction  between  the  components 


of  the  random  process  at  X.  ^ and  X , then  the  value  of  bicoherence  is 


equal  to  one  and  the  biphase  is  the  phase  lag  of  the  nonlinear  process, 


where  <f> ( X ^. ) are  the  phases  of  the  interacting  components  of  the  random 


process . The  case  in  which  «^(X  , X ) = 0 or  iris  called  phase  locking 


by  Barnett  et  al.  (19  71  ). 


2.  Definitions  of  Cross-bispectrum  Density,  Cross -biphase,  and 
Cross- bi  coherence 


The  definitions  concerning  cross-bispectrum  have  been  derived 


following  the  same  line  of  the  bispectrum  analysis  by  Hasselmann  ( 1966) 


and  Roden  and  Bendiner  (1973).  For  stationary  random  functions  of 


time  (or  distance)  with  zero  mean,  sayx(t),  y(t),  z(t),  the  cross- 


ll 


1 


spectrum.  F (\  ) , and  cross-bispectrum  densities,  B^iX  , ■ t-  2>) 
are  defined  respectively  as  the  Fourier  transforms  of  R (T).  co- 
variance  of  x (t ) and  y(t);  and  S cross-bicovariance  of 


x(t),  y(t)  and  z(t): 


F (X  ) 
xy 


r 

2 IT  J 


-iZirXT  , 

R (T)e  dT 

xy 


where 


R (T)  = < x(t+T)y(t)>  * 
xy 


B (X  , X ) = 
xyz  1 2 


r\  00  - 2 n(X  T.  + ^2^2  ^ 

-Tiff  V 

(2ir)  -oo 


where 


s (T,.  T,)  = < X(t+T.  )y(t+T?)z(t)> 

XV  z 1 2 i 2 


In  terms  of  the  components  of  the  Fourier  Stietjes  representa- 


tion of 


x<ti ff 


dZ  (X  )e 
x 


i2irX  t 


y (t  ) = j dZy(X  )e 
- oc 


i2irXt 


z(t)  = \ dZ^(X  )e 


i2trX  t 


< dZx(X1)dZ  (X2)  > = Fxy(X  i >dX  J . if  Xj+X2  = o 


= 0 ifXj  + X2  0 


12 


i.  e.  <dZ  (X  )dZ  *(\  ) > = F (Y  )d\ 
xi  y 1 xy  1 1 


< dZx(X  , )dZy,X  2)dZ(X  3)>  = B^ft  J , X 2)dK  [dX  2 it  X 1 + X 2n  3 =0 

= 0 if  \ j + \ 2 + \ ;£  0 (II 


< dZx (X , )dZy(x 2,dZz *,x ; + X 2)> - Bxyz (Xj.X^dXj 


dX. 


The  symmetry  relations  are 


R (T)  = R (-T) 
xy  yx 


Sxyz(Tl’  T2)  Syxz(T2’  T1 ] ~ Syzx(T2'Tl  ’ *T1  } 


' Szyx(-T2’  VV 


and  hence 


VM  = V‘M 


Bxyz(XrK2)  Byxz(X2’>ll)  Byzx(X  2’  'X  1~X  21 
= Bzyx(*XrX2’Xl>  * 


(24 


(25) 


26) 


( 27) 


The  physical  meaning  of  cross- spectrum  density  is  well  known. 
It  is  a complex  quantity.  Its  amplitude  shows  whether  a frequency 
component  in  one  random  process  (say  x(t))  is  associated  with  larger 
or  smaller  amplitudes  at  the  same  frequency  in  other  random  process 
(say  y (t)).  The  phase  shows  whether  the  frequency  component  in  one 


random  process  lags  or  leads  the  component  at  the  same  frequency  in 
the  other  random  process.  As  pointed  out  by  Roden  and  Bendiner  (1973), 
the  cross-bicorrelation  measures  the  joint  interaction  between  three 
sets  of  variables  displaced  relative  to  each  other  in  the  time  domain 
and  gives  an  indication  of  the  persistence  of  the  interaction.  The 
cross-bispectrum  density  is  a measure  of  such  interaction  in  frequency 
domain.  If  any  two  of  the  three  random  processes  are  the  same,  then 
the  cross-bispectrum  density  gives  the  measure  of  the  quadratic  inter- 
action between  two  frequency  components,  say  \ , \ in  one  random 
process,  and  the  frequency  component  of  - (\ X ^ ) in  the  other  ran- 
dom process.  If  all  three  are  the  same  process,  the  cross-bispectrum 
density  is  just  the  bispectrum  density  of  the  process;  subsequently,  all 
the  subscripts  in  the  formula  can  be  dropped. 

Cross-bispectrum  density  is  also  a complex  function.  Its  phase 
can  be  defined  the  same  way  as  the  biphase: 


<(>  (X.  . . \ = tan 

xyz  1 2 


1 ImlBxvz(W) 


RetB*yhW] 


(28, 


where  4>xyz.(X  j , \ is  called  cross -bi phase  by  Roden  and  Bendiner 


(1973). 


By  the  same  analogy,  the  cross -bicoherence  can  be  defined  as 

2 f 20) 


->  B (L.lJdLdL, 

bic2  ft, . 1 z — '—Jj 

xyz  1 2 _ 3 

(Fx(VFy(X2)Fz(VdX  > 
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_2££Ti 


' 


where  b»cxy>z  j . ^2)  is  the  cross -bicoherence,  B (X.^,  \2)  is  the 
cross -bispectrumdensity  and  FJ\  j )•  Fy  (X  £ j,  FJX  J are  the  spectrum 
densities  of  random  processes  x(t),  y(t)  and  z(t)  at  frequencies  \ . 

X2’  x 3 respectively.  The  cross-bicoherence  is,  by  definition,  a 
real  quantity  equal  to  or  larger  than  zero.  In  cases  where  the 
nonlinear  process  is  due  to  quadratic  interaction  between  the  indepen- 
dent frequencies  X } and  X ^ such  that  the  component  at  X.  ^ + X is  equal 
to  or  proportional  to  the  product  of  the  components  at  X ^ and  \ , the 
cross -bicoherence  (and  also  the  bicoherence)  is  unity.  It  follows, 
therefore,  that  if  the  cross-bicoherence  differs  from  unity  signifi- 
cantly, the  nonlinear  processes  cannot  be  of  the  simple  type  (Roden 
and  Bendiner,  1973). 

The  interpretation  of  cross-biphase  can  also  be  obtained  by  the 
analogy  to  that  of  biphase.  In  general,  one  can  consider  that 
oispectrum  analysis  is  the  special  case  of  cross -bispectrum  analysis 

B-  Derivation  of  Rotary  Bispectrum  and  Rotary  Cross-bisoectrum 
F unctions 

1.  Fourier  Representation  ofRotary  Components  of  a Two-dimensional 
Vector  Random  Process 

The  idea  of  rotary  components  was  introduced  by  Mooers  (1970). 
His  method  is  based  on  a decomposition  of  the  complex- valued  series 
into  polarized  components.  For  a two-dimensional  velocity  vector 
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expressed  as  a complex  valued  quantity,  oscillatory  elements  execute 
periodic,  elliptical  orbits  in  the  hodograph  plane.  The  trajectory  of 
the  velocity  vector  in  this  plane  then  is  composed  of  a superposition 
of  those  ellipses  of  all  oscillatory  elements.  The  velocity  vector 
series  can  be  decomposed,  for  each  frequency,  into  two  counter- 
rotating circular  components.  Each  component  has  its  own  amplitude 
and  phase.  The  counter-clockwise  rotating  components  correspond  to 
oscillatory  motions  at  positive  frequencies;  the  clockwise  rotating 
components  correspond  to  those  at  negative  frequencies.  In  other 
words,  the  oscillatory  element  at  each  frequency  is  composed  of  two 
counter-rotating  components  which  have  angular  velocities  equal  to 
the  positive  and  negative  angular  frequencies  respectively  (Gonella, 
1972). 

The  oscillatory  elements  can  be  represented  as  the  Fourier 
components  of  the  velocity  vector  series  (Perkins,  1970).  Gonella 
(1972)  gives  the  Fourier  transform,  for  the  angular  velocity  uj  of  the 
complex  time  series  u(t)  as 


T 

U ( w)  = — (*  u(t)e  lwtdt 

r\ 


U (u)  = U(w)  e 


where  T is  the  duration  of  the  complex  time  series  u(t),  and 


L6 


u(t)  = Uj  (t)  + iu2(t) 


uj  and  u2  are  the  scalar  components  of  u1  along  horizontal  rectangular 
axes  0Xj,0x2.  The  complex  coefficient  U(w)  gives  the  amplitude 
j U(w)  | and  phase  at  the  initial  time,  of  the  rotary  components 
with  an  angular  velocity  w.  The  angular  velocity  w will  be  an  integral 
multiple,  positive  or  negative,  of  2tt/T.  U(w)  can  be  readily  obtained 
in  terms  of  the  sine  and  cosine  Fourier  coefficients  of  u^  (t)  and  u2(t) 
corresponding  to  the  angular  frequency  a = | co  | , 


1 fT 

7 J [ ui  (t)  + iu2(t)  ] [cosut  - i sin  wt  ] dt 


[Aj(«)  + B2(w)  ] + i[A2(u)  - (w)  ] 


where 


i rT 

A. (u))=—  \ u.  (t)coscotdt 

J J0  J 

i rT 

B.  (w)  - — \ u.  (t)sinwtdt 

.1  J o J 

A and  B are  the  sine  and  cosine  Fourier  coefficients  of  u (t). 
J J j 


U +(<r  ) - U (u>)  if  id  = a 
U (ff  ) = U (w)  if  id  = -a 


Then 


U+(ff  ) = [Aj  (a  ) + B2(ff  )]  + i [ A 2 (cr  ) - B1  (ff  )] 
u.(ff  ) = [ Aj  (<r  ) - B2((r  ) ] + i [A 2 (ff  ) + B}  (cr  )] 
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U + (o-  ),  U (cr  ) are  complex  quantities.  They  represent  the  amplitude 
and  the  phase  of  the  rotary  components  of  angular  frequencies  +a  and 
-a  respectively.  They  are  called  rotary  Fourier  coefficients  in  this 
study.  It  is  obvious  that  the  rotary  Fourier  coefficient  of  (+0-  ) does 
not  necessarily  equal  the  conjugate  of  the  rotary  Fourier  coefficient 
of  (-c  ) . 

The  inverse  of  Fourier  transform  will  be 


u(t)  = u1  (t)  + iu^(t)  = 2 


+ 00 


U (w)e 


i cot 


(33) 


n = - 00 


where 


(d  = 


2irn 


In  terms  of  rotary  Fourier  coefficient 


v-»  00 

Uj  (t)  + iu^ (t)  = ^ U+(ff  )e1<7  + U (o'  )e  'l  °”  4 . (34) 

^ o-=0 

If  U1  o'  U2o-  ^ represent  the  scalar  components  of  the  oscillatory 
element  with  angular  frequency  o-  , then  the  ellipse  equation  is 


Iff 


(t) 


+ iu2ff  (t)  = 


U+(ff  )e  ^ ^ 4 + 


U (ff  )e 


-iff  t 


(35) 


By  equation  <32)> 


2 

one  can  see 


The  amplitude  of  rotary  component  derived  in  this  study  is  twice  of 
that  given  by  Gonella  (197  2). 
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(3 


uj  a (t)  = 2[  Aj  (a  )cosiTt  + B,  (cr  )sina  t J 
U2a  ^ ~ A2^<7  ^COS  17  1 + ^s'ncr  1 ] 

The  rotary  Fourier  .-('efficients  derived  here  are  the  linear 
functions  of  the  sine  and  cosine  coefficients  of  the  scalar  components 
of  the  vector  time  series.  Since  the  scalar  components  are  resolved 
along  an  arbitrary  chosen  -ectanglar  axes,  their  Fourier  coefficients 
are  not  independent  of  the  coordinates.  In  fact,  if  one  rotates  the 
coordinate,  the  Fourier  coefficients  will  ooey  the  same  coordinate 
transform  as  the  scalar  components  (Fofonoff,  1969),  i.e.: 

A,  = Aj  cos  6 + A 2 sin  0 

A ' = -A.  sin  6 + A_  cos  0 

L-  L £ 

(37 

' = B cos  0 + B?  sin  0 

B2'  = -B  sin  0 + B,  cos  0 

wnere  0 is  the  angle  (measured  counterclockwise)  between  the 
primed  and  unprimed  axes.  Substituting  equation  (37)  into  equation 

( 3 2)  yields 

U+'(ff  ) =U  +(cr  )e"19 


From  equation  '),  it  is  obvious  that  the  amplitude  of  the  rotary 


Fourier  coefficient  is  invariant  under  the  coordinate  rotation  but  the 
phase  does  change  by  an  amount  equal  to  the  angle  of  the  rotation; 


U±'(a  ) = U±(cr  > 


= $ - 9 

± (r  ± cr 


In  the  terminology  of  an  ellipse,  U (cr  )|  + U (a  ) is  the  length 


of  semimajor  axis  of  the  ellipse,  and  U+(o-  ) - U (o-  ) is  the 


length  of  its  semiminor  axis.  The  orientation  of  the  major  axes, 
from  the  data  coordinate  system,  is 

$ . $ 

% 


and  that  of  the  minor  axes  is 


$>  + <j> 

ff+  * 


(Wooers,  1970).  Therefore,  if  |u+(<r)|  equals  | U (cr ) I the  oscillatory 
element  is  in  rectilinear  oscillation;  if  either  u+(°’  ) or  (J  (cr  )| 
is  zero,  the  oscillatory  element  describes  pure  circular  motion 
rotating  in  the  direction  of  the  non- zero  component.  For  a scalar 
random  process,  the  rotary  Fourier  coefficients  will  be  the  same  as 
the  ordinary  Fourier  coefficients  (see  equation  32).  In  this  special 


U+(ff  ) = U_  U ) . 


case,  and  the  oscillatory  element  is  also  a rectilinear  oscillation. 


2.  Advantages  of  Using  Rotary  Component  Method  for  Spectrum 
Analysis  of  Two-d  imens ional  Vector  Random  Processes 


f 

i' 


! 


For  a two-dimensional  vector  random  process,  spectrum  analy- 
sis can  be  conducted  on  the  scalar  component  series  resolved  along 
rectangular  coordinates.  The  spectral  functions  generated  from  the 
two  scalar  series  are  dependent  upon  the  orientation  of  the  coordinates 
along  which  the  observations  were  made  or  analyzed.  Also,  the  two 
component  series  are  not  independent.  These  problems  were  noted 
by  Fofonoff  (1969).  He  proposes  that,  for  each  frequency,  a normal 
set  of  coordinates  can  be  introduced  such  that  the  non-invariant 
cospectrum  between  the  scalar  component  series  vanishes.  This  is 
a tedious  method  for  analyzing  one  set  of  such  vector  random  process 
alone,  to  say  nothing  of  analyzing  a pair  of  them.  The  rotary  com- 
ponent method  can  overcome  these  problems  as  shown  by  Perkins 
(1970),  Mooers  (1970,  197  3)  and  Gonella  (1 972). 

Since  certain  oscillatory  components  in  the  ocean  current  are 
expected  to  have  regular  motion  with  one  sense  of  rotation,  for 
instance  inertial  oscillation  and  tidal  oscillations,  this  method 
appears  more  adequate  to  analyze  the  oceanic  current  velocity  data 
than  the  resolved  scalar  component  series  method.  The  rotary 
spectral  functions  and  rotary  cros s - s pectral  functions  have  been 
well  defined  and  discussed  by  Gonella  (1972)  and  Mooers  (197  3). 

Some  of  these  functions  will  be  reviewed  briefly  where  the  rotary 


i 

■ 


i 

JL 


: 
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component  method  is  extended  to  third  order  spectral  analysis. 


3.  Rotary  Spectrum-  and  Rotary  Bispectrum-Functions 

In  earlier  sections,  the  spectrum  ana  Dispectrum  densities  have 
been  defined  as  the  frequency  decomposition  of  a real  random  process' 
second  and  third  moments  respectively.  For  a vector  random  process, 
the  same  operations  yield  general  second  and  third-order  spectrum  den 
sities  (Brillinger  andRosenblatt,  1966b).  If  the  decomposition  is  carried 
out  in  terms  of  angular  velocities  on  the  covariance  and  bi covariance 
functions  of  a two-dimensional  random  process,  the  results  are 
defined  as  the  rotary  spectrum  density  and  the  rotary  bispectrum 
density  respectively.  Let'u(t)  be  the  stationary  two-dimensional 
vector  random  process.  Its  covariance  is  defined  as 

R(T)  = < u (t)u(t  + T)  > 

. 

= < u (t)u  (t  - T)  > (41 ) 

(Bendat  and  Piersol,  1966).  By  the  same  token,  its  bicovariance  is 
defined  as 

SfTi » T2)  = < u*(t)u(t  + Tj  )u(t  + T2)  > 

= < u (t)u(t  + T-,)u(t  + Tj  ) > . (42) 

The  second  and  third-order  spectrum  densities  of  u (t)  are  then 
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P(w)=j^y  R(T)e'lwTdT 


j pp”  -*  "i(wlTl +w2T2) 

RBluru2):~IJ]  S(We  dTidT2 

\CJ1)  -OC 


where  w is  the  angular  velocity  of  an  oscillatory  component. 
Expressed  in  terms  of  the  Fourier-Stieltjes  representations, 


*•>-2 


tt/  \ •,.3  2,-irn 

U (w)e  , with  w = ~ ^ — 


< U (u>j ) U (w^)  > = P(w)dw  if 

= 0 if  w ^ w 


* 2 

< U (u>j ) U (w^)  U (u>3)  > = RB(uy  w^dw  if  ^ +(*>2  = u» 


if  Wj  tui^  i 


For  u>  = ± cr  , equation  (45)  is  rewritten  as 


Pu±u±  (<r  )d<1  = <U±(OU±  (»)>*<  U ± (ff  ) > (4V) 


where  P (a  ) and  P (a  ) are  the  rotary  spectrum  densities  of 
u+u+  u-u-  if 


In  practice,  w i9  limited  by  the  Nyquist  frequency,  Ny,  i.  e.  , 

- 1 n_  1 TT  TT 

2At  - T - 2At  ’ "At-W-At 

where  At  is  the  data  interval. 


MM 


u(t)  for  positive  and  negative  angular  velocities  respectively.  In  view 
of  the  characteristics  of  the  rotary  Fourier  coefficients  given  in  the 
previous  section,  it  is  obvious  that: 

1)  u are  ^dependent  of  the  coordinates  along  which 
the  scaler  components  of  u(t)  are  measured;  and 

2)  .,,,(0)  is  not  s'inai  to  P (ct)  in  general,  except  when  the 

U+U+  u-u-  r 

angular  frequency  w is  zero  or  the  amplitudes  of  the  rotating 
components  are  equal  (rectilinear  oscillation);  and 

3)  the  sum  of  P (cr)d(j  andP  (c)da  is  the  contribution  to  the 

u+u+  u-u- 

total  variance  of  u(t)  from  the  oscillatory  element  of  the  frequency4  ] a| . 

In  case  of  the  third -order  spectrum  density,  all  three  angular 
velocities  can  have  both  positive  and  negative  values.  The  constaint 
relationships  of  uk  , j = 1,  2,  3,  have  eight  permutations: 


Gonella  (1972)  gives  a complete  set  of  functions  of  the  rotary  com  - 


ponent method  for  spectrum  and  cross -spectrum  analysis.  P (o)do 

u±  ui 


defined  here  is  twice  the  value  of  his  mean  kinetic  energy  spectrum 
S±  (a ): 


S±(  a)  = 1/2  < U±(ct)U±  *(CT)  >. 

He  also  shows  that  all  the  rotary  spectral-  and  rotary  cross  spectral 
functions  can  be  related  linearly  to  the  ordinary  spectral-  and  cross 
spectral  functions  of  the  resolved  scalar  component  series. 

It  should  also  be  noted  that  common  usage  in  spectrum  analysis 
is  to  refer  to  the  second -order  spectrum  estimates  as  "energy"  even 
though  the  described  quantity  does  not  have  "energy"  dimensions. 

This  conventional  terminology  is  also  adopted  in  this  study  e.g.  , 
the  energy  spectrum  of  wind  stress. 
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K 


3>-’l  + u 2 = ^ 3 

4)  <r  1 

5)  Ci  -a2  = -3 

6 ) -<T 

7> 

8)  -<r 

-r.,^3 


2 = _<r  3 


1 - a 2 = 3 


(48) 


The  sum  of  these  permutations  can  be  seen  clearer  in  a frequency 
domain  diagram  of  Euclidean  representation  /Fig.  1). 

Due  to  the  symmetric  relation  of  ^ and  (Equation  46')  there 
are  only  four  independent  permutations,  i.e.,  the  third  spectrum  den- 
sity of  a two-dimensional  vector  random  process  can  be  defined  in 
the  half-plane  of  die  frequency  domain.5  Correspondingly,  the  rotary 
bispectrum  densities  covering  the  four  independent  sections  in  the 
frequency  domain  are 


The  third  order  spectrum  density  of  the  two  dimensional  vector  ran- 
dom processes  can  also  be  represented  by  the  linear  combination  of 
ordinary  bispectrum  density  and  cross-bispectrum  density  of  the  vector's 
resolved  scalar  component  series.  There  are  eight  independent  bispectrum 
and  cross-bispectrum  densities  for  each  such  combination.  It  is  not 
only  very  tedious  to  compute  them  but  also  very  difficult  to  interpret 
them. 
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RB^, cruder  =<  U + (<x  { )U  + (a2  )U+*(<r  3>  > , a +a  = 


RB(-cr1,a2)dT  = < UJdl)U  + (cr2)U  + "'t(<r3)>>  -aj +a2  = 


cr  _ 

> a . 

2 

1 

(49) 


RB^j.  '^2)dcr  =<  U+(<rl)U_(cr2)U_JfC(cr3)>>  a } - cr  = <r  . 


RB('a]’  '<r2,da  = < u _ (cr j )U  (cr 2 )TJ _ ^(a- 3 ) > , -<r  1 -cr  = -a 


Since  the  rotary  bispectral  density  is  a complex  quantity,  it 
can  be  expressed  in  polar  form,  for  instance: 


R B (<r  j , cr  2 ) = R B (cr  j , cr  2 ) e 


iR^fo-j , o-2) 


(50) 


where  R4>(<r^,a  ) is  defined  as  rotary  biphase: 

(m[  RB(aJ , a 2)J 


R<Mo",  , o-  ) = tan 

12  Re  [ RB(aJ , cr2)]  _ 


(H  l) 


The  corresponding  bicoherence  is  defined  as  the  rotary  bicoherence, 


Rbic(<x  ^ , a 2 ): 


Rbic  (o-j  , cr2)  = 


RB(tr  j , a 2)do-  2 


21  2 


(Pu+u+(ai,Pu+u+(<r2)Pu+u  + (<r3,dff  > 


(5  2) 


It  is  easy  to  see  that  the  amplitude  of  the  rotary  bispectrum  density 
and  the  rotary  bicoherence  are  independent  of  the  rectangular 
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- 


i 3 


coordinates  along  which  the  scalar  components  of  the  vector  random 


process  are  measured.  However,  the  rotary  biphase  is  not  an  inde- 
pendent quantity.  If  the  coordinate  system  is  rotated  by  some  angle, 
the  value  of  the  rotary  biphase  will  change  the  same  amount  in  the 
direction  opposite  to  that  of  the  coordinate  rotation.  This  can  be 
verified  by  using  equations  (38),  (49),  (50)  and  (51). 

The  rotary  bispectrum  is  then  the  measure  of  low  order  non- 
linear interactions  among  the  rotary  components  of  the  oscillatory 
elements. 


4.  Rotary  Cros  s - spectrum  Functions  and  Rotary  Cross-bispectrum 
Functions 


The  rotary  cros  s - spectrum  and  rotary  cross-bispectrum 
functions  are  obtained  here  following  the  same  line  of  reasoning  as  in 
the  previous  section.  Letx(t),  y(t)  and  "z(t)  be  the  stationary  two- 
dimensional  random  processes.  The  cross -covariance  is  defined  as 
(Bendat  and  Piersol,  1966). 


R (T)  = < x*(t)y(t+T)>  = < x(t)y*(t-T)  > <53) 

xy 


and  the  cross -bicovariance  is  then 


“S  (T.  , T,)  =<  x (t)y(t+T,  )z(t+T,)> 
xy  z 1 <L  1 c 

= <x* (t)z*(t+T_  )y(t+Tj ) > 


= S (Tv  T.  ) 
xzy  2 1 


» oo 


P~<“>-£  I «jT)e-‘^dT 


xy 


- oo 


xy 


-»  ■ —z  N 1_.fr..  T,).-"1'*  ”2'2 


" X w,  » I 

xyz  1 2 


(55) 


(2ir) 


- oo 


00  - -i(w,T,  +w,T,) 

xyz'V  '!'''  dTldT2  (56) 


where  PxyM  is  the  cross  - spectrum  density  between  x(t)  and  y(t),  and 


xy 

"xy z "°1  ’ ‘*'2' 

Expressed  in  terms  of  the  Fourier-Stieltjes  representations 


’ w">^  18  t*ie  cr°s s - bi spectrum  density  among  x(t),  y(t)  and  z(t). 


x(t) 


■l 


oo 


n = -oc 


U (w)e 
x 


iu>t 


y(t)  = y U (co)elwt 
^n=-oo  y 


z(t) 


■I 


00 


n - - oo 


u Me 

Z 


icjt 


< Ux(wJ,Uy  (w2}  >_  PxyS,dw  if  W1  = w, 


= 0 


if  W1  * w2 


(57) 


(58) 


<Ux(“l)Uy("2,Uz  =RBxyz(“l'"2)d“2  « ",  +“2  - u, 

if  Wj  +w2  t (59-) 


= 0 


it  is  obvious  that, 


and 


p M = p M 

xy  yx ' ' 


(W1  > <*>_)  - RB  (w  , w ) • 
xyz  l i xzy  2 1 


(60) 
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For  lo  = ± u , equation  (55 ) is  rewritten  as 


pxiyl('»d'  = <U*±(,|V  M> 


(61) 


where  Px+  +(«r ) and  Px_y>)  are  the  rotary  cross  spectrum  densities 
between  x(I)  and  y(t)  for  the  positive  and  negative  angular  frequencies 

respectively. 

In  the  case  of  rotary  cros  s-bispectrum  density,  the  permutation 
and  symmetric  relations  are  the  same  as  those  for  rotary  bispectrum 
density.  There. 

RBxyz(<ri,(r2)dCr  2 =<  Ux+(<T1  )Uy+(Cr2)Uz+  ^ > 

* 1 +ff2  =a3 


1 


RBxyz('Crrtr2,d<T  2 " < Ux-^l)Uy  + (tr2)Uz+  (<T3):> 


RBxyz^l’  _cr2)d<J  =<Ux  + (al)Uy-(<r2)Uz-  (ff3)- 


ffrff2  = _<r3  a2 


> o' 


(62) 
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and 


RB  (-<r  -«r  )da  2 = <UX.(^1)U  >2>Uz-  (<r  3f  > 

Xyz  16  7 


i - 2 = -ff  3 


RBxyz(ffl,<T2>  = Ky^r0^!6 


iR.4>  (cr  , . <r  J 

xyz  1 £■ 


(63) 
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where  (“"j  > o^)  is  de*ined  as  rotary  cross -biphase 

Im[  RB  (a  ,j  )] 
R*  - tan'1  *rz,  1 2- 

xyz  1 2 RetRB™z 'V'z'] 


Correspondingly,  the  rotary  cross-bi coherence  is 


Rbic  (<t  , o-  ) = 
xyz  1 2 


RB  (<r  ,<r  )da 
xyz  1 2 


(Px+x  + (VPy+y+<VPz+z^3)d<r\ 


Again,  one  can  see  that  the  rotary  bispectrum  functions  are  the 
special  cases  of  the  rotary  cross-bispectrum  functions.  If  all  thr< 
vector  random  processes  are  identical,  all  the  subscripts,  x,  y,  z, 
can  be  dropped,  and  one  has  a set  of  rotary  bispectrum  functions. 
But  the  rotary  cross-biphase  has  a higher  dependence  on  the 
coordinates.  It  is  easy  to  see 


R'4>  = R<t>  - 6 - 6 + G 

xyz  xyz  xyz 


where  9 , 0 , 0 are  the  amount  of  rotation  of  which  the  coordinates 
y z 

ofx(t),  y(t),  z(t)  undertake  respectively. 

The  rotary  cros s -bispectrum  is  the  measure  of  low  order  non- 
linear interactions  among  the  rotary  components  of  the  different 
vector  random  processes.  If  the  first  two  processes  are  identical, 
then  the  rotary  cross  bispectrum  is  the  measure  of  quadratic  inter- 
actions between  the  two  different  processes. 
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5.  Estimation  Procedure 

The  rotary  bispectrum  and  rotary  cross-bispectrum  functions  can 
be  easily  estimated  with  the  rotary  Fourier  coefficients  of  the 
respective  two-dimensional  vector  random  processes.  The  algorithm 
is  readily  supplied  by  the  related  equations  given  in  the  previous  sec- 
tions. The  proposed  procedures  are  outlined  and  given  in  the  Appendix. 

6.  Applicability 

The  rotary  bispectrum  is  a reliable  tool  for  examining  the  non- 
Gaussian  nature  of  a vector  random  process.  Being  independent  of 
the  rectangular  coordinates  along  which  the  scalar  components  of  the 
process  are  described,  the  amplitude  of  the  rotary  bispectral  density 
gives  the  unique  measure  of  the  low-order  nonlinear  interactions 

I 

(the  quadratic  interactions)  internal  to  the  process  itself.  However, 
the  amplitude  does  not  show  the  cause  and  effect  of  the  quadratic 
interaction  as  it  is  only  a measure  of  the  contribution  to  third  moment 
of  the  process  by  combinations  of  three  oscillatory  interacting  com- 
ponents whose  frequencies  are  related  in  the  manner  described  by 
equation  (48),  i.e.,  one  frequency  is  the  algebraic  sum  of  the  other 
two  frequencies.  It  does  not  show  explicitly  the  amount  of  internal 
energy  transferred  among  the  oscillatory  components  (e.g.,  the 
third  moment  of  a current  series  has  dimension  of  L3/T3,  which  is 

-3 

I 


32 


not  an  energy  dimension).  The  rotary  biphase  shows  the  phase  lag 
of  the  quadratic  interaction.  Since  the  rotary  biphase  is  not  a unique 
quantity,  i e. , not  independent  of  the  coordinates,  it  does  not  yield 
meaningful  information  if  considered  alone. 

The  amplitude  of  the  rotary  cross  - bispectrum  density  a quantity 
also  independent  of  coordinate  description,  gives  a unique  measure  of 
the  quadratic  interactions  among  three  different  vector  processes  at 
the  triple  frequencies  of  which  one  is  the  algebraic  sum  of  the  other 
two.  If  a pair  of  the  three  processes  are  the  same  process,  then  the 
measure  is  of  quadratic  interaction  between  the  two  independent 
processes  (in  this  case,  two  of  the  three  oscillatory  elements  belong 
to  one  process).  Since  the  rotary  cross -bispectrum  analysis  involves 
two  or  three  different  vector  process,  the  role  of  the  cause  and 
effect  of  the  quadratic  interactions  may  be  assigned  to  each  process 
by  examining  the  physical  nature  of  the  interacting  process.  One 
assigns  the  'cause"  to  the  two  oscillatory  components  and  the  "effect" 
to  that  component  whose  frequency  is  the  algebraic  sum  of  the  fre- 
quencies of  the  'cause"  components.  The  possible  combinations  of 
such  an  assignment  are  illustrated  in  Fig.  2. 

In  this  study,  the  quadratic  interactions  to  be  examined  are  the 
effect  of  two  interacting  oscillating  wind  stress  components  upon  the 
energy  contained  in  an  ocean  current  at  the  third  frequency.  Thus, 
one  can  postulate  the  relation 


I 
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XlLLs. 


a)  [±ffI]+[±^2]=[±^3]  b)  [±^2]  + [±ff3]  = [±^]  c)  [±«r1]  + [±ff3]  = f±a2] 

0<a,  < N 0 < <r  _ < N 0 < ct  „ < N 

1-y  — 2 — y ~ 3 - y 

Figure  2.  Possible  combinations  of  cause  and  effect  for  quadratic 
interactions  shown  by  the  rotary  cross -bispectrum. 


Path  1 Path  2 

Figure  3.  Quadratic  energy  transfer  paths  from  process  Is  to 
process  c. 


c(t)  = as(t)  + hs  (t) 


(67) 


! 

I 

f 
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where  s(t)  and  c(t)  are  the  wind  stress  and  current  fields  respectively, 
a is  the  coefficient  of  the  linear  dependence  term  and  b is  the  coeffi- 
cient of  the  quadratic,  nonlinear  dependence  of  current  on  the  wind 
stress.  Two  cases  are  examined  in  a subsequent  technical  report; 

1)  wind  stress  vs.  current  at  14m,  and  2)  wind  stress  vs.  current 
at  34  m.  In  a similar  manner,  quadratic  interactions  between  the 
ocean  currents  at  14  m and  34  m levels  are  examined  (in  this  case, 
the  left  hand  side  of  the  equation  (67)  will  denote  current  field  at 
34  m and  the  right  hand  side  will  denote  current  field  at  14  m). 

Thus,  the  rotary  cross-bispectrum  analysis  involves  only  two  pro- 
cesses of  which,  in  case  of  wind  stress  vs.  current,  the  wind 
stress  field  is  assigned  the  role  of  "cause",  and,  in  the  case  of 
current  vs.  current,  the  current  field  at  14  m is  assigned  that  of 
"cause"  process. ^ 

The  amount  of  energy  transferred  still  cannot  be  shown  quanti- 
tatively as  the  rotary  cross-bispectrum  density  only  shows  the  amount 
of  common  cube  shared  by  the  three  interacting  components,  and  its 
dimension  is  the  product  of  the  energy  of  the  "cause"  process  and  the 


Different  method  for  examining  the  nonlinear  relationship  is  given  by 
Dewan  (1969).  However,  his  method  only  deals  with  the  interaction 
among  the  two  components  with  same  frequency  in  one  process  and 
third  component  in  another  process,  whose  frequency  is  either  zero 
or  twice  of  that  of  the  other  two  components.  Thus,  it  only  covers 
the  axes  of  a 3 = 0,  a ) - a 2 in  the  frequency  domain  shown  by  Fig.  1 . 
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half  power  of  the  "effect"  process.  Besides,  the  energy  transformed 
quadratically  from  one  process  to  another  can  take  two  paths:  1) 
direct  quadratic  interaction  among  the  two  components  in  the  "cause" 
process  and  the  interacted  component  in  the  "effect"  process,  2) 
linear  interaction  between  an  "effect"  component  and  a "cause"  com- 
ponent which  is  itself  produced  by  quadratic  interaction  internal  to 
the  "cause"  process.  These  two  paths  are  shown  in  Fig.  3.  They 
may  both  exist  in  nature  but  the  cross -bispectrum  cannot  show  the 
latter  even  qualitatively,  in  order  to  assess  the  amount  of  energy 
transformed  quadratically  from  one  process  to  another,  an  additional 
scheme  is  developed  in  the  following  section. 


C.  Energy  Transfer  Functions  between  a Pair  of  Two-dimensional 


Vector  Random  Processes 


1.  Rotary  Fourier  Representation  of  Quadratic  Interaction  between 
Two  Different  Vector  Random  Processes 


Hasselmann  (1966)  developed  a Fourier  representation  of  the 
quadratic  interactions  between  two  random  processes.  A modifica- 
tion makes  it  useful  for  the  rotary  component  representation  of  two 
different  vector  random  processes.  The  rotary  Fourier  coefficients 
previously  defined  are  used  in  place  of  the  Fourier  Stieltjes  represen 
tations.  Also,  where  Hasselmann  used  the  ordinary  conjugate  to 


denote  a function  of  negative  frequency,  the  rotary  component  method 


: > 

:•  t 


I < 
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requires  use  of  the  function  of  negative  frequency  itself.  With  this 
modification,  Hasselmann's  equation  (7)  is  applied  to  the  rotary 
ponent  analysis,  as 


com- 


V'< 


icr  t 


[uc  + ((T)elu  1 + U Jo-  )e~ia  *] 


((T)Us+(«T)ela  1 +T_(a)U  _(a)e'1<rt] 


'II  + ,68) 


al  a 2 


K>  i * 2,Us-(<T1  )us.^2)e"i(<Tl+<T2,t 


+ V(<r  -o  )U  (<r  )U  (o  , Je*  1~  ° 2 }t 
1 £ s t i s - C 


V(-2«r1)Ue.„i,Us+«r2,ei(-<rlW2''i 


where 


1 1 Us±  (cr ) , lrc±  (a  ) are  the  rotary  Fourier  coefficients  of  random 
processes  s(t),  u (t ) respectively,  e.g.  wind  stress  and  current 
velocity . 

2)  a is  the  frequency  with  a range  of  0 < a < x . For  the  dis- 
crete time  series,  the  upper  limit  of  a is  the  Nyquist  frequency. 

3»  T±(cj)  are  the  time-independent  linear  transfer  functions, 

Sjc 

an<3  T±  (o  ) ^ T (o  ) in  general. 
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4)  K (ff  +<r  ) are  the  time-independent  quadratic  sum  transfer 

± 1 2. 

functions  and 

K±(<V*2)  = K±(Vrl) 

but  K±(<r1+a2)  * K±*(a1+(T2)  in  general. 

5)  V(«r  -»r2)  and  are  the  time  independent  quadratic 

difference  transfer  functions. 


V(«r1-ff2)  = V + (ff1-ff2) 


if  «r1-ff2>  0 


= V (a  1 - cr 2 ) if  o-j-a  2<  0 
V(-aI+a2)  = V + (-a1-«r2)  if  ' V°  2 > ° 

= V (-<J1+ff2)  if  -cr  1 +cr  2 < 0 

Thus  it  is  obvious 

v±(a1-a2)  = V±(-a24v1) 

and  V+(<r1-<r2)  / V ^a1-a2)  in  general. 

For  any  particular  frequency  a 3 equation  (38)  can  be  separated 


into  two  different  equations  concerning  positive  angular  velocity  ^ = 
a and  negative  angular  velocity  ^ -a  jComponent.  • Thus  if  | *2  |> 


[ f 


+ 2 > > V (-a  +a  )U  (<r  )U  (a  »ei(-°l+a2,t 
* “ s ~ 1 s + 

x 6[  (‘<J1  +a2]  ' a 3 j C69) 


ffl  a2 


and 


Uc-(<73)e  Cr3t  = T_K)U  - (<r,)e'1Cr3t 


3'  s-v  y 


+ 2^  y K 


ZZ 

<rj  tr2 

-i  (cr  +cr  )t  / 
e 1 2 6 [ - (cr  j +a  ] 


(70) 


+ 2 


n 


V (0-  -a)  U (a  ) U (<r:,)e'l(<rrlI2,t 
1 s+i  s-2 


a (j 


6 [-(«r1-ff2)  + *,] 


where  6[  ] is  the  Dirac's  Delta  function  serving  as  a constraint 

device  for  the  double  summation  e.g.,  S^f^)  - a ] implies  the 
summations  only  cover  the  range  . Because  of  this 

straint,  these  two  equations  can  be  simplified  as 


con- 


Uc*<°3>  stS> 


"ll 

J1  °2 


K*<V‘’2)Us*(VUs±,‘r2, 

x6[±(a1  + o-2)±(7  ] 


(71) 
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V1  W",±<V  usi(»2) 


°1  a2 


x 6 [ ± (a  l - a 2 )±  a ] 
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where 


One  has  to  notice  also  that  each  U (<r  ) may  be  expressed  in 

s±  3 

equation  (71)  form  if  there  are  any  quadratic  interactions  existing  in 

u (t).  However,  for  simplicity,  it  will  be  assumed  that  u (t)  is  a 
s s 

Gaussian  process.  But  this  assumption  shall  be  corrected  later  when 
the  quadratic  transform  functions  K and  V are  defined  (see  Hassel- 
mann,  1966). 

2.  The  Transfer  Functions 


The  transfer  functions  T,  K,  V,  can  be  defined  in  terms  of 
the  Rotary  spectrum  parameters. 

a.  Linear  Transfer  Functions  T (cr ) 


The  linear  transfer  function  is  given  by  Hasselmann  (1966)  as 
the  ratio  of  the  cross  - spectrum  density  between  the  driving  process 
and  the  responding  process,  to  the  spectrum  density  of  the  driving 
process.  For  rotary  spectrum  expression,  one  has 


^Dimensional  wise.  [T  (<r)l  =[U  (a)U  (<r)~*l 

tK*<VV)  =|V±<*»,*o2)]  = [U  c±(.3)U>±(.l)-‘  U („  I'1] 


(72) 


P («r) 
T (a)  =~-S-+C+- 
+ Ps+s+(ff) 


and 


T fu)  = 


P (<r) 

s-c- 

P (*> 

s - s - 


(73) 


b.  Quadratic  Sum  Transfer  Functions,  K ( cr  + cr  ) 

±1  2' 


Multiplying  the  positive  frequency  equation  of  equation  (71)  by 

* * 

U s+(ffj  ) and  U 3-^2^’  w^t^<rj+0’^=ff2’  anc^  taking  the  ensemble 


average,  one  has 
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RB 


ssc(Vff2)d(j2  = T + (ff3,RB5f!Sss(<7l+<T2)dcr2 

+ 2K  + (a1hr2)Ps+s+(a1)Ps+s+(a2)da 


(74) 


If  cr  is  zero  frequency,  a 3 = a q = 0.  Then 


^^ssc^O’  2 = W’^sssW'*"2 

+ 4K  + („o^2)psts  + (»o)Pst5 +(02)d»' 


3 

In  a manner  analogous  to  that  given  in  Hasselmann  (1966), 

< U (cj.  ) U (w  )U  *(W,)U  '"(w. ) > = 0 

si  s 2 s 3 s 4 

unless  Wj  = ijj u>2  = . If,  for  the  rotary  component  method 


(75) 


°2  = ^4  = ~ a 2’  wi  ~ w3  =<7i’  then 


< US4,VUs>2,OVUs-’(<V>  = IWVPs.,.('’2>d‘' 


as  K + (<r0+<72)  = V+(-cr  and  Pg+s+(^0)  = P8_s>0)-  Here  the 

* 

rotary  bispectral  density  RB  (<r,,c r_)  is  included  to  account  for 

sss  1 2 

the  non-Gaussian  nature  of  the  u(t)  process  itself.  Solving  for  the 
quadratic  sum  transfer  function,  one  obtains 


RB  (cr  ,a  ) 

v , , v ssc  1 2 

K+  ff  1 + 2 = T~P 


(ff,)P  (cr  _ ) 
s +s  + 1 s +s  + 2 


V*3> 


RB  (°^  , cr  ) 
sss  1 Z 


3'2Ps+s+((VPs+s+<‘V* 


+ O'  . 


= a 


(76) 


The  first  and  second  terms  on  the  right  hand  side  of  equation  (7^) 

_y 

account  for  the  path  1 and  path  2 illustrated  in  Fig.  3.  If  u(t)  is 

indeed  a Gaussian  process,  then  the  second  term  does  not  exist 
* 

because  RB  (cr  , cr  ) will  not  be  significantly  greater  than  zero. 

On  the  other  hand,  if  there  are  no  quadratic  interactions  between  the 
two  random  processes,  but  the  driving  process  has  its  own  internal 
quadratic  interaction,  then  the  cr  component  of  the  responding 
process  still  will  be  affected  by  cr , cr^  components  in  the  driving 
process  provided  that  the  linear  transfer  function  T+(cr^)  is  signifi- 
cant. The  quadratic  sum  transfer  function  is  then 

* 


K + (ffl  + °2>  = ~ 


This  interaction  is  not  explicit, 
negative  components. 


(<J3)  2P 


RB  (cr.  ,cr  ) 
sss  1 2 


s+st'VWV 


The  same  relation  also  holds  for 
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K ( a + or  ) 


RB  (-cr  , -o'  ) 

ssc  1 2 


1 2 2P  (cr  )P  (a  ) 

s-s-  1 s-  s - 2 


- T (<r) 


RB  (-cr  , -o-  ) 
. sss  1 2 


3 2P  (a,  )P  (o-  ) 
s-s-  1 S-S-  2' 


(77) 


c.  Quadratic  Difference  Transfer  Functions,  V (t  cr  ± it  ) 

± 1 2 


The  quadratic  difference  transfer  functions  can  be  derived 
the  same  way  as  the  quadratic  sum  transfer  functions.  Thus 


RB'  (=Fcr  , ±o_) 

V fro  ±o  ) = i L. 

* 1 2 2P  (ff,  )P  (o  ) 

"""  ~ T 1 S ± s ± 2 


STS 


RB  (To  , ±o  ) 

, V sss  1 2 

Wlp ~ 75 7T7  • (78) 


3 2P  (o  )P  (o  ) 
st st  1 s±  s±  ' 2 


3.  The  Energy  Transfer  Equati 


on 


imen- 


The  energy  transfer  equation  between  a pair  of  two-di 
sional  vector  random  processes  can  be  readily  obtained  from  equati 

'3'  ' c±  '"3'  ~c±  '“3 

2 


on 


7711.  As  pc±c±(a Jdo  = < U (a  ) U ___  (oo)>,  one  has 


Pc±c±(ff3)dff  = 


T (o  ) 
± 3 ' 


P (o  Jdo  + 

s±  s±  3 


Mathematically,  there  should  be  a factor  of  four  in  the  second  and 
thira  terms  on  the  right  hand  side  of  equation  (78)  to  account  for  the 
permutations  of  each  pair  of  the  interacting  frequency  components. 
The  factor  of  four  is  not  included  in  this  equation  as  only  one  of  the 
permutations  represents  the  true  energy  process. 


t d 
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K*(<rl+ff2,|2  Ps±s±(<ri)PS±s±(iI2)dcr2  ’ 6 l 1 ±0l±(r2)=F(T3  ] 


v±  (=F^±a2) 


2 PsTST(al,PS±S±((r2)d(r2  *6[  <«1*^2)»»3] 


(79) 


al  az 

where  ja^|  > |cr  | . 

The  first  term  on  the  right  hand  side  of  the  equation  denotes  the 
amount  of  energy  being  transformed  from  the  driving  process  to  the 
responding  process  through  linear  action.  The  second  and  third 
terms  denote  the  amount  of  energy  being  transformed  through  quad- 
ratic interactions. 


D.  Statistical  Considerations 


The  development  of  this  proposed  technique  is  based  on  two 
hypotheses;  first  that  the  random  processes  are  not  normally  dis- 
tributed (non-Gaussian),  and  second  that  the  random  processes  are 
stationary.  Therefore,  any  set  of  data  used  for  this  analysis  should 
first  pass  a statistical  test  of  these  hypotheses.  Since  this  analysis 
is  carried  out  on  a limited  sample,  significance  test  should  also  be 
applied  to  function  estimates  such  as  rotary  cross-spectrum  density, 
rotary  bispectrum  density  and  rotary  cross-bispectrum  density.  The 
estimates  should  be  significantly  greater  than  zero  to  yield  any  use- 
ful information.  These  tests  are  discussed  in  the  following  sections. 
For  detailed  testing  procedures  the  reader  is  referred  to  referenced 


literature. 


1 . Normality  Test 


The  well  known  X test  for  goodness  of  fit  of  normality  is 
chosen  as  the  main  testing  scheme.  It  has  been  well  documented 
(Crow,  Davis  and  Maxfield,  I960),  and  the  proper  computer  pro- 
grams needed  for  the  test  are  also  available  (Bevington,  1969).  A 
complimentary  testing  scheme  using  Cornu  ratio  and  skewness 
(Crew  and  Bodvarsson,  1971  ) will  also  be  used  for  comparison. 

The  tests  are  applied  to  the  scalar  component  series  for  the  two 
dimensional  vector  random  process.  If  both  of  the  scalar  component 
series  fail  all  of  the  tests  then  this  vector  random  process  is  not 
normally  distributed  (Jenkins  and  Watts,  1968). 

2.  Stationarity  Test 

For  practical  reasons,  the  rigorous  testing  of  stationarity  of 
a random  process  is  not  feasible  (Bendat  and  Piersol,  1966).  How- 
ever, in  this  proposed  analysis,  all  of  the  data  series  used  must 
satisfy  the  condition  of  weak  stationarity,  to  second  order,  as  a 
minimum  condition.  This  is  due  to  the  fact  that  all  the  functions 
defined  here  are  derived  directly  or  indirectly  from  the  Fourier 
coefficients  of  the  scalar  component  series.  It  is  rather  important 
to  ensure  that  the  Fourier  representation  of  the  series  is  valid. 

Haubrich  (1965)  gives  such  a test  in  detail,  based  on  the  idea 
that  the  coherence  squared  between  the  Fourier  coefficients  at 
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different  frequencies  will  be  significantly  greater  than  zero  if  the 
random  process  is  not  stationary.  This  idea  has  been  further  proven 
by  Lumley  (1970).  He  states  that,  for  a stationary  process,  Fourier 
components  at  two  different  frequencies  are  uncorrelated,  though  not 
statistically  independent. 

Therefore,  the  stationarity  test  by  Haubrich  (1965)  is  utilized 
for  this  analysis  for  the  sake  of  compatibility  and  practicality.  It 
will  be  conducted  on  each  resolved  scalar  data  series  to  ensure  that 
its  Fourier  representation  is  valid. 

% 

3.  Significance  Test  of  the  Rotary  Cross -spectrum  Estimates 


The  rotary  cross -spectrum  density  is  meaningful  only  when  its 
corresponding  coherence  estimate  is  significantly  greater  than  zero. 


The  rotary  coherence  y^ficr)  *s  defined  as  (Wooers,  1 97  3) 


y\z  P 


P1±2±"> 


1±1±(<7)  P2±  2±  (<T) 


(80) 


where 


o<y12  <1  • 


and  Pj±|±(°')  is  the  rotary  cross- spectrum  density  between  processes 
1 and  2,  and  P (o  ),  P (<t  ) are  rotary  spectrum  densities  of 
processes  1 and  2 at  ±cr  . It  provides  a non-dimensional  measure 
of  the  correlation  between  two  time  series  as  a function  of  frequency. 
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Jenkins  and  Watts  (1968)  give  a simple  significance  test  scheme 


for  such  correlation  based  on  least  squares  analysis  in  the  frequency 
domain.  They  show  that  the  random  variable 


(EDF-2)  Y122(<0 

2(1  - \12V)) 


(81 ) 


is  distributed  approximately  as  the  F where  EDF  is  the 

[JLjLJJO  ■ l j 

equivalent  degree  of  freedom  of  the  spectrum  estimates  F 

2,  (EDF-2) 

is  the  Fisher  F distribution  with  2 and  (EDF-2)  degrees  of  freedom, 
andY^  is  t^e  estimated  value  of  the  coherence  squared.  Thus, 


(EDF-2)y122  (ff) 

2(1  -?  *(„»  £ (EDF-2)'1 


(82) 


then  y is  not  greater  than  zero  at  the  a significance  level. 

This  test  is  derived  from  ordinary  spectrum  and  cross- spectrum 
functions.  It  can  also  be  applied  to  the  rotary  component  spectrum 
functions  by  following  the  same  analogy. 


4.  Significance  Test  of  the  Rotary  Bispectrum  and  Rotary  Cross-bispectrum 
Estimates 


The  third  order  spectrum  estimates  of  a random  process  (vector- 


valued or  real -valued)  are  asymptotically  unbiased  and  asymptotically 


normal  under  mild  conditions  (Brillinger  and  Rosenblatt,  1966a).  Thus 

2 

the  estimates  of  the  coherence  squared  asymptoti cally  has  a X distri- 
bution with  two  degrees  of  freedom  since  it  is  sum  of  squares  of  two 
terms  each  of  which  is  asymptotically  normally  distributed. 

Haubrich  (1965)  gives  that  the  expected  values  of  the  estimate  of 
the  bicoherence  squared,  obtained  from  truncated  data,  is  approxi- 
mately 2/EDF  if  its  corresponding  bispectrum  density  estimates  have 
true  vanishing  value. 

Therefore,  the  bispectrum  density  estimates  will  be  significantly 
greater  than  zero  at  the  5%  level  if  its  correspinding  estimates  of  the 
bi coherence  squared  is  greater  than  or  equal  to  6/EDF.  This  testing 
criterion  can  also  be  applied  to  the  cros s -bispectrum  estimates  which 
have  been  shown  as  the  general  case  of  the  bispectrum  estimates  in 
earlier  section.  By  the  same  analogy,  it  can  also  be  applied  to  the 
rotary  bispectrum  and  cross-bispectrum  estimates.  Thus,  if  the  esti- 
mate of  the  rotary  bicoherence  squared  is  less  than  6/EDF,  the 
corresponding  rotary  bispectral  density  is  not  significantly  greater 
than  zero  at  the  5%  level.  The  same  rule  holds  for  rotary  cross- 
bispectum  density. 

Due  to  the  fact  that  the  bicoherence  can  only  be  estimated  when 
the  three  interacting  components  all  have  a reasonable  amount  of 
spectrum  density  (Hinich  and  Clay,  196 8),  a significance  test  should  be 
applied  to  the  spectrum  density  estimates  of  the  rotary  components 
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before  their  bispectrum  density  are  tested  for  significance.  One  can 
see  that  this  requirement  is  also  true  for  the  significance  test  of  the 
rotary  cross-bispectrum  density  estimates. 

5.  Significance  Test  of  the  Rotary  Spectrum  Density  Estimates 

For  a discrete  time  series,  the  spectrum  density  estimates 
should  decrease  to  zero  before  the  Nyquist  frequency.  If  the  spectrum 
density  estimate  at  the  Nyquist  frequency  appears  greater  than  zero, 
this  estimated  value  is  the  noise  introduced  by  aliasing  of  the  fre- 
quency components  above  the  Nyquist  frequency.  This  noise  cannot  be 
completely  eliminated  in  practice.  Heuristically,  the  noise  is  re- 
garded as  the  systematical  error  for  the  true  zero  spectrum  density 
estimate,  and  it  can  be  used  as  a threshold  value  for  testing  other 
non- zero  spectrum  estimates.  The  confidence  limits  for  the  spectrum 
density  estimate  at  a certain  significance  level  can  be  easily  obtained 
with  the  estimate,  its  equivalent  degrees  of  freedom,  and  the  variable 
associated  with  the  significance  level  (Jenkins  andWatts,  1968).  A sim- 
ple test  is  then  proposed  in  this  study  to  compare  the  lower  limit  of 
any  rotary  spectrum  density  estimate  against  the  upper  limit  of  that 
at  the  Nyquist  frequency.  If  the  former  value  is  less  than  the  latter, 
the  tested  rotary  spectrum  density  is  considered  not  to  be  significantly 
greater  than  zero  at  that  certain  level. 

It  is  obvious  that  all  of  the  significance  tests  are  essential  to 


the  reliability  of  this  proposed  analyzing  technique.  Thus,  in  order 
to  have  valid  results,  all  the  interpretations  and  calculations  should 
only  be  based  on  the  rotary  spectrum  parameters  which  have  passed 
these  tests. 

This  process  has  been  applied  to  a set  of  wind  and  current  data 
gathered  from  the  TOTEM  buoy  in  1970.  The  data,  the  analyses,  and 
the  results  are  given  in  Bispectrum  and  Cross -Bispectrum  Analysis  of 
Wind  and  Currents  off  the  Oregon  Coast:  II  Data  and  Examples. 
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Estimation  Procedures  of  the  Rotary  Spectral,  Rotary  Cross- 
spectral,  Rotary  Bispectral  and  Rotary  Cross-bispectral 

Functions 

A . Editing  the  Data 

For  practical  reasons,  the  existing  data  of  two-dimensional 
vector  random  process  are  often  recorded  in  the  form  of  a finite 
length  discrete  time  series  of  their  resolved  scalar  components  or 
of  their  magnitude  and  directions.  Since  rotary  Fourier  coefficients 
are  obtained  from  the  linear  combinations  of  the  cosine  and  sine 
Fourier  coefficients  of  the  scalar  components  and  since  their 
amplitudes  are  independent  of  the  coordinates  along  which  the  scalar 
components  are  resolved,  data  should  be  in  the  form  of  resolved 
scalar  component  series  along  E-W  (ox^)  and  N-S  (ox^)  directions 
with  E,  N designated  as  the  positive  axes.  The  choice  of  coordinates 
is  arbitrary;  to  be  consistent,  conventional  notation  is  used.  For 
the  sake  of  saving  computation  effort,  the  cosine  and  sine  Fourier 
coefficients  are  computed  with  the  help  of  the  relatively  new,  fast 
Fourier  transform  algorithm  (Welch,  1967).  Also,  in  order  to 
increase  the  number  of  degrees  of  freedom  of  the  estimates  of  various 
functions,  (i.e.,  to  ensure  the  proper  statistical  stability),  each  data 
record  is  broken  into  several  truncations,  from  which  the  averaged 
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estimated  functions  are  obtained  (Haubrich,  1965;  Hinich  and  Clay, 
1968;  Huber  et  al. , 1971).  Therefore,  the  data  records  have  to  be 
edited  to  achieve  the  proper  estimation  with  less  computational 
expense. 

1 . Rearranging  the  Data  Interval  and  Removal  of  the  Mean 

The  data  sampling  interval  should  be  of  proper  length  so  that 
the  Nyquist  frequency  will  not  be  excessively  higher  than  the  upper 
limit  of  the  frequency  range  of  interest.  This  can  be  done  with  either 
a block  average  of  several  closely  spaced  data  values  or  with  a 
weighted  moving  average  technique.  The  latter  method  is  a combi- 
nation of  low  pass  filtering  and  averaging.  With  properly  chosen 
weighting  coefficients,  the  aliases  can  be  avoided.  However,  certain 
data  points  on  both  ends  of  the  record  will  be  sacrificed.  Then  the 
mean  of  the  data  should  be  removed  in  order  to  be  consistent  with 
the  assumption  that  the  random  process  has  zero  mean. 

2.  Adjusting  the  Data  Record 

The  data  record  may  have  to  be  adjusted  to  provide  the  proper 
number  of  truncations  so  that  the  averaged  function  estimates  can 
have  proper  equivalent  degrees  of  freedom.  Because  of  the  require- 
ments of  the  Fast  Fourier  transform,  the  number  of  data  points  in 
each  truncation  must  be  a power  of  two.  In  order  to  have  the  proper 


number  of  truncations,  it  may  be  necessary  to  adjust  the  data  length 
by  either  adding  zeros  at  the  end  of  each  truncation  (Henry  and 


Graefe,  1971)  or  by  overlapping  a part  of  data  in  each  truncation 
(Haubrich,  1965;  Welch,  1967;  Huber  etal,,  1971),  The  latter 
method  is  proposed  in  this  study  because  the  equivalent  degrees  of 
freedom  EDF  of  the  spectral  estimates  from  the  overlapping  trun- 
cated data  can  be  calculated  explicitly  (Welch,  1967).  Let  x(t), 
t = 0,  ...  N-l  be  the  finite  discrete  time  series  of  a random  process 
x(t).  If  x(t)  has  a record  length  of  T and  a data  interval  of  At,  the 
total  number  of  data  points  N is  T/At.  Overlapping  can  be  done  by 
taking  truncations  of  length  mAt  (m  data  points  in  each  truncation) 
with  the  starting  point  D data  points  apart.  Then  the  truncated  data 
is  represented  as 

x(k,  j ) = x(j  + (K-  1)D)  j = 0,  ...  N-l 

k = 1 , . . . K (A  - 1 ) 

where  k is  the  total  number  of  the  truncations. 

Ks<2-»i  + i (A  -2> 

There  are  M-D  data  points  overlapped  in  each  truncation.  If  no 
overlapping  is  allowed,  one  simply  lets  D equal  M.  Then  equation 
(A  • 2)  becomes  K = N/M  (Welch,  1967). 
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3.  Applying  weighting  function  to  the  truncated  data 

This  modification  can  reduce  the  leakage  of  spectral  power 
from  a spectral  peak  to  frequencies  far  away  (Haubrich,  1965; 
Welch,  1967;  Huber  £t  ad. , 1971).  However,  weighting  the  data  not 
only  will  degrade  the  frequency  resolution  but  also  will  change  the 
expectations  and  the  variances  of  the  spectral  estimates  (Huber 
et  al.  , 1971).  Since  the  data  available  for  this  study  are  already 
limited  in  length  this  technique  will  not  be  used  for  the  sake  of 
keeping  adequate  frequency  resolution. 

B.  Approximating  the  Rotary  Fourier  Coefficients  of  the  Truncated 
Data 

For  each  truncation  of  M data  values,  (1  + M/2)  number  of  the 
rotary  Fourier  coefficients,  with  the  frequencies  ranging  from  zero 
up  to  Nyquist,  can  be  obtained  from  the  linear  combinations  of  the 
sine  and  cosine  Fourier  coefficients  of  the  resolved  scalar  compon- 
ents in  that  truncation. 

Let  u (k,  j)  and  u (k, j)  be  the  scalar  components  of  the  kth 
1 £ 

truncation.  Then  the  discrete  fast  Fourier  transform  is 

z<k,X)^y  u (k,j)W(j)  e'l2lTjq/M 

P ML  - = Q P 

M 

p = l,2  ; q = 0 . . . . , — (A  -3) 

(Welch,  1967) 
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where,  Zp(k,  X)  is  the  complex  Fourier  coefficient;  \ is  the  frequency 

of  cycle  per  data  interval  CPH  and  X.  = q/M;  W (j ) is  the  weighting 

function  or  data  window.  If  no  weighting  function  is  used,  then 

W (j  ) = 1 , for  j = 0,  ....  M- 1 . Physically  (k,  X.  ) represents  the 

amplitude  and  the  phase  of  the  qth  frequency  component  of  the  scalar 

components;  the  resolved  frequency  band  of  the  component  is  1 /MAt 

10 

CPD  if  no  data  window  is  applied  . The  cosine  and  sine  coefficients 
of  Up(k,j)  are  the  real  and  imaginary  parts  of  Z^(k,  X.)  respectively 
by  definition.  Thus 


A (k,  X ) = Ref  Z (k,\  )1 
P 1 P J 


Bp(k,  X ) = Im[  Zp(k,  X)] 


(A  -4) 


By  equation  (32),  the  rotary  Fourier  coefficients  of  the  truncated 
data  u(k,  j ) are 


U+(k,  X ) = [ Aj  (k,  X ) + B2(k,  X)]  + i [ A (k,  X ) - B (k,  X )] 


U _ (k,  X ) = [Aj  (k,  X ) - B2(k,  X )]  + i [ A2(k,  X ) + B (k,  X )] 


(A  -5) 


The  definition  of  the  discrete  Fourier  transform  is  not  uniform  in 
the  literature^.  Some  authors  use  a quantity  equivalent  to  Z(X  )MAt 
or  Z(X  )(MAt)  , others  use  a positive  exponent  instead  of  a negative 
one  (Jenkins  and  Watts,  1968;  Hinich  and  Clay,  1968).  These  variants 
appear  not  only  in  magnitude  but  also  in  physical  meaning.  For 
example,  the  quantity  Zp(X  )(MAt)  represents  the  distribution  of  up(t) 
per  resolved  frequency  band,  and  Zp(X  )(MAt)1/<2  represents  the  dis- 
tribution of  up(t)  per  square  root  of  the  resolved  frequency  band. 

For  the  purpose  of  consistency  with  the  derivation,  the  coefficients 
of  DFT  defined  as  given  in  equation  (A -3). 


C.  Computing  the  Spectrum  Estimates 


Once  the  rotary  Fourier  coefficients  are  calculated,  computing 
the  spectrum  estimates  is  fairly  straightforward.  Using  the  definitions 
given  before,  the  formulae  for  computing  various  spectrum  estimates 
are  given  in  the  following  equations. 

For  a single  vector  random  process,  x(t) 

,K 
'k=l 


P (X)d\ 
x±x± 


— y1 

K L, 


U (k,  X) 


(A  -6) 


RB(±\  , ±X  2)dX  = 
K 


K 4=1 


Ux±(k'S)  Ux±(k'K2)Ux± 


±X1±X2  = ±S 


(A  -7) 


RB(±\  f ±X  2)  = 


rb(±x1>±\2) 


eiR4>(X1,  \2)  (A  _g) 


Rbic  (±X  , ±X2> 


* 

2 

K 

l “xl»'S>U,i<kA2>U*t  <k’S» 

^k=l 

r 

^k=i 


2VK 

2 

ux±(k,  Xj) 

Ux±(k.X2) 

u (k,X  ) 

x±  3 

(A  -9) 


where  A denotes  estimating  value 
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11More  sophisticated  average  procedures  are  suggested  by  Huber  et 
al.  (1971  ) for  their  study  of  ordinary  bispectrum  analysis.  Bust  these 
procedures  need  too  lengthy  a data  record  and  too  much  computing 
time  to  be  used  in  this  study. 
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In  the  case  of  multiple  random  processes,  let  U (k,  X ),  U (k,  X ), 

U (k,  X ) represent  the  rotary  Fourier  coefficients  of  their  truncated 
z± 

data.  Then 


P . ,(k)dX 
x±y± 


U Ik,  X ) U . (k,  X ) (A  -10) 

x±  y± 


RB  (±  X , ± X )dX 
xyz  1 2 


kL  , Ux±(S)Uy±(XZ)Uz± 

k=l 


tXl±X2=±S 


(A  -11) 


Rbic  ±\  , ± \ ) 
xyz  1 i. 


i 


jj. 

K Z Ux±(k’Xl)Uy±(k,X2)Uz± 


2f 

u .(k,  X ) 

v±  c 

2 ^ 
l 

-K 

> 

U , (k,  X ) 
z±  3 

^k=l 

7 

■"k  = l 

(A  -12) 
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